% This script calculates the return a prospect theory investor is willing
% to forgo using NYSE/AMEX/Nasdaw returns.
clear all;

addpath('../Matlab/');
format long;

% prospect theory parameters
lambda = 1;
delta = 1;
alpha = 0.5;

% other parameters
window = 1; % size of windows of returns to look at
skipoct87 = 0;
readdata = 1;
bins = 10; % break data into deciles
plotfigure = 0; % plot figures

if readdata == 1
    % read in daily returns
    fileID = fopen('daily returns NYSE.csv','r');
    keepreading = 1; % when not end of file

    row = 1;
    valid_ids = [];
    while keepreading
        % read in a line
        line = textscan(fileID, '%f %f',1,'HeaderLines',1,'Delimiter',',');
        % check if we've reached the end of the file
        if isempty(line{1})
            keepreading = 0;
        else
            % extract the values we need, converting to numeric
            date(row) = double(line{1});
            daily(row) = double(line{2});
        end

        row = row+1;
    end
    fclose(fileID);
    save('calibration');
else
    load('calibration','date','daily');
end

if skipoct87 == 1
    date = date(date < 19871001 | date > 19871100);
    daily = daily(date < 19871001 | date > 19871100);
end

% identify monthly boundaries
monthidx(1) = 1;
idx = 2;
currentmonth = floor(mod(date(1)/10000,1)*100); % extract month
for i = 2:size(daily,2)
    % identify start of month
    month = floor(mod(date(i)/10000,1)*100); % extract month
    if (month ~= currentmonth)
        monthidx(idx) = i;
        idx = idx + 1;
        currentmonth = month;
    end
end
nummonths = size(monthidx,2);
monthidx(nummonths+1) = size(date,2)+1; % add index for end of data 

% get 6 month returns
midx = 1;
for i=1:nummonths-5-window % need 6 months for return plus window of data after that
%     if mod(i,6) == 1 || mod(i,6) == 7
%         sixmonth(midx,:) = [i prod(daily(monthidx(i):(monthidx(i+6)-1))+1)]; % gross return
%         midx = midx + 1;
%     end;
    sixmonth(i,:) = [i prod(daily(monthidx(i):(monthidx(i+6)-1))+1)]; % gross return
end

% get window of daily returns following 6 month periods with lowest and highest
% returns
numsixmonths = size(sixmonth,1);
numinquantile = floor(numsixmonths/bins);
sixmonth = sortrows(sixmonth,2); % increasing in return
postlow = [];
posthigh = [];
for i=1:numinquantile
    postlow = [postlow daily(monthidx(sixmonth(i,1)+6):monthidx(sixmonth(i,1)+6+window)-1)];
    posthigh = [posthigh daily(monthidx(sixmonth(numsixmonths-numinquantile+i,1)+6):monthidx(sixmonth(numsixmonths-numinquantile+i,1)+6+window)-1)];
end

for j=1:3
    if j == 1 % all data
        disp('All data');
        data = daily;    
        [allmean, allstd, allskew, allmaxneg, allindiff] = bootstrap(data,delta,lambda,alpha)
    elseif j == 2 % post low returns
        disp('After low returns');
        data = postlow;
        [lowmean, lowstd, lowskew, lowmaxneg, lowindiff] = bootstrap(data,delta,lambda,alpha)
    else
        disp('After high returns');
        data = posthigh;
        [highmean, highstd, highskew, highmaxneg, highindiff] = bootstrap(data,delta,lambda,alpha)
    end
    
    if plotfigure ==1
        % plot EU as a function of private signal on mean return
        returnrange=-0.01:0.0001:0.01;
        n = size(data,2);
        idx = 1;
        for i=1:size(returnrange,2)
            buyEU(idx) = EU_la(1/n*ones(n,1),data'+returnrange(i),delta,alpha,lambda);
            sellEU(idx) = EU_la(1/n*ones(n,1),-data'-returnrange(i),delta,alpha,lambda);
            idx = idx + 1;
        end
        figure;
        plot(returnrange,buyEU,returnrange,sellEU);
        legend('Buy EU','Sell EU');
    end
end
    